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ABSTRACT: We propose a model for diffusion-limited annihilation of two species, A + 
B — > A or B, where the motion of the particles is subject to a drift. For equal initial 
concentrations of the two species, the density follows a power-law decay for large times. 
However, the decay exponent varies continuously as a function of the probability of which 
particle, the hopping one or the target, survives in the reaction. These results suggest that 
diffusion-limited reactions subject to drift do not fall into a limited number of universality 
classes. 
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Diffusion-limited reactions (DLR) in low dimensions have been extensively studied be- 
cause of their anomalous kinetics, and because they serve as simple prototypes of complex, 
nonequilibrium dynamical systems [1-10]. Recently, attention has focused on anisotropic 
DLR: the particles move preferentially in one spatial direction, giving rise to an effective 
drift [11-16]. 

For the one-species DLR of coalescence; A + A — > A, and annihilation; kA — > 0, it 
has been shown that the spatial anisotropy introduces no appreciable changes [11, 15-16]. 
However, the two-species annihilation process, A + B — > 0, exhibits markedly different 
behavior with or without drift [12-14]. In one dimension, when the initial concentrations 
of A and B particles are equal, the total concentration decays as t -1 / 4 in the isotropic 
case, and as t -1 / 3 in the presence of drift. This interesting behavior is attributed largely 
to the hard core repulsion between particles of like species. A heuristic explanation, based 
on the Burgers equation, for the t -1 / 3 decay has been presented along these lines [13-14]. 
A recent exact analysis of a two-species annihilation model with no hard core interactions 
confirms their essential role: in their absence, the concentration decays as t -1 / 4 , with or 
without the drift [16]. 

In this work, we introduce a new anisotropic two-species annihilation model, with 
hard core interactions between like species. A distinct characteristic of our model is that 
it contains an adjustable parameter. We find that the concentration power-law decay 
exponent depends continuously on this parameter. Nonuniversal exponents are a surprising 
feature. Our finding is based on extensive numerical Monte Carlo simulations to be detailed 
below. 

Our model is defined on a one-dimensional lattice. Each lattice site can be in one 
of the three states: empty (0), occupied by a single A particle, or occupied by a single 
B particle. Numerical results were obtained for the case of maximum anisotropy [12-14], 
where particles may hop only to the nearest site to their right: 

■■■AO >---0A---, ■■■BO >---0B---, (1) 
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with equal rates for A and B. Hopping is disallowed if the target site is occupied by 
a particle of the same species; this models hard core interactions. If the target site is 
occupied by a particle of the opposite species, hopping is allowed and reaction takes place 
with outcome determined by the probabilistic rule: 

••••»'••• •{::::;»::: :::!::? „-•••»•»••• •{:::;;';::: „ ® 

The reaction may then be represented as 

A + B ^ Aor B, (3) 

where the product species is that of the hopping particle, with probability p, or that of 
the target particle, with probability 1 — p. Thus, p may be thought of as a parameter 
which represents the "persistence" of the hopping particle. Stoichiometrically, however, 
the reaction rule is equivalent to the symmetric annihilation reaction 

A + B^±A+±B, (4) 

regardless of the value of p, because the number of AB nearest-neighbor pairs is equal 
to the number of BA nearest pairs (assuming homogeneous initial distributions), even for 
unequal concentrations of A and B particles. 

We have performed concurrent simulations of the above model on a cluster of over 
50 IBM RS6000 workstations, for lattices of 10 6 sites, with periodic boundary conditions. 
Simulations ran up to times t = 10 6 . Particle hopping was random and independent, 
with the rate defined so that in one time unit each of the surviving particles performs 
an average of one hopping attempt. The attempt could lead to motion, reaction, or be 
discarded, depending on the state of the target site. Data from about 50 such runs were 
collected for five different persistence values: p = 0, 0.25, 0.5, 0.75, and 1. All runs 
began with an 80% full lattice, with equal concentrations of A and B particles, distributed 
randomly. 
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In Figure 1, we show a log- log plot of the concentration decay as a function of time 
as obtained from simulations. The data exhibit an asymptotic power-law decay with a in- 
dependent exponent. In Figure 2, we plot the local slopes of the decay curves as a function 
of time. The decay exponents are obtained through extrapolation to time t — > oo. Our 
estimates for the decay exponent for different values of the persistence parameter are listed 
in Table 1. 

One of the recognized fluctuation effects in two-species reactions is the tendency of 
the system to develop large alternating domains of A and B particles [7,13-14,17-18]. We 
also collected data for some quantities of interest associated with this phenomenon. These 
include: (^ab(£)), the average distance between nearest particles of opposite species; 
(lAA(t)), the average distance between nearest particles of the same species; (L(t)), the 
average length of a domain of like particles, measured between the first particles of adjacent 
domains; and (n(t)), the average number of particles per domain. In the large-time 
asymptotic limit, each of these quantities grows as a power of t. The relation c ~ (n) / (L) 
was used to check for consistency. Simulation results (from about 600 runs on 10 5 -site 
lattices) are summarized in Table 1. The error margins in the table are meant only as 
rough estimates and may be overly optimistic. They merely represent the spread of values 
we have obtained by extrapolating the data with different powers of 1/t. 

The case of p = 1 seems special in that the number of particles per domain remains 
bounded and of order unity, even after very long times. In addition, the concentration 
decays as t -1 / 2 , similar to isotropic one-species coalescence, A + A — > A. These results are 
particularly interesting in view of the following considerations. 

Suppose that through fluctuations, large domains of A and B particles were created, 
for general p. Left alone, each domain moves with the drift velocity, while the particle dis- 
tribution in it is governed by the hard-core constraint and is related to the Burgers-equation 
theory; see [13-14] and references therein. However, reactions at both ends modify the do- 
main structure. Our results suggest that the interplay of the two processes is nontrivial 
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for all < p < 1, and the domain size measures involve nontrivial, varying exponents; see 
Table 1. Self-organization by formation of large domains has been viewed as an essential 
condition for non-mean-field behavior in low-dimensional reactions [6-7,17-18]. 

However, for p = 1, we find that there are no large domains. This indicates that 
if a large domain of, for instance, A particles, were formed, it would be destroyed by B 
particles "catching up" from the left. Indeed, for p = 1 the incoming, hopping particle 
determines the species of the reaction product with probability 1. Thus, regardless of the 
size of the A-domain, its denser, left end will be "eaten up" by B particles faster than its 
dilute, right end could catch up with other A particles, while eliminating the B particles 
in between. The formation of the shock-like profile due to the interplay of hard-core and 
anisotropy is essential for this phenomenon, i.e., we need the two domain ends to have 
different densities. 

Once the p — 1 system is mixed, one would naively anticipate the mean-field behavior, 
~ 1/t. However, we know from reactions like A + A — > 0, or A, that non-mean- field 
fluctuations can also arise in "well-mixed" situations (no domains), in the form of non- 
mean-field interparticle distributions [19]. Detailed studies of the interparticle distribution 
for the p = 1 model will be reported in the future. Since all known one-dimensional 
reaction-diffusion systems, with or without hopping isotropy, possess only one "diffusive" 
length-scale, ~ t 1 ^ 2 , the exponent values of 1/2 for the p = 1 model (see Table 1) are 
probably exact. 

The case of p = 1/2 resembles the anisotropic two-species annihilation model, A+B — > 
0, studied by Janowsky [12-13] and by Ispolatov et al. [14]. Notice that the difference in 
the number of particles between the two species, Na — Nb, is locally conserved (in every 
reaction event) in the Janowsky [12] model, but only globally conserved in the present 
model. This distinction is of no serious consequence in the isotropic case of no drift [20] . In 
the case of drift, we argue that the two models should be similar. In particular, the heuristic 
arguments of Ispolatov et al. [14] could be repeated for our model with p = 1/2, almost 
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without change, obtaining the same conclusions and predictions. Indeed, our simulation 
results for p = 1/2 are not inconsistent with those already published for A + B — > 0. 

An interesting observation is that for p = 1/2 the convergence to the long time asymp- 
totic behavior of our model seems faster than that of the Janowsky model. An example is 
shown in Figure 3, where we plot the time-dependent exponent of the concentration decay 
(obtained from local slopes of log(c) vs. log(t)), for both cases. We have no explanation for 
this phenomenon. Perhaps it could be exploited to settle the slight discrepancies between 
the numerical findings of Janowsky [12-13] and the results of Ispolatov et al. [14]. 

In summary, we have proposed a one-dimensional DLR model with continuously vary- 
ing exponents. It is important to emphasize that our density-exponent values were deter- 
mined largely by numerical studies, and furthermore, the actual estimates range from 1/2 
for p — 1 down to near 1/3 at p = 1/2, and to about 1/4 for p = 0. Since all three values: 
1/4, 1/3, 1/2, and only these values, have been encountered for various universality classes 
of other one-dimensional reactions, one could suspect that our observation of continuous 
exponents is an artefact of finite-time numerics. However, our numerical simulations were 
really "large-scale" by modern standards, and the data presented seem to suggest error- 
limits which clearly favor continuous variation rather than a smooth crossover between 
three universal values. 
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CAPTIONS 



Figure 1: Concentration decay for different values of the persistence parameter. The solid 
curves represent simulation data for p = 0, 0.25, 0.5, 0.75, and 1 (top to bottom). 

Figure 2: Local power-law decay exponents for p = 0, 0.25, 0.5, 0.75, and 1 (top to 
bottom) as a function of inverse time. The power-law decay exponent are obtained 
from the local slopes of the curves in Figure 1. 

Figure 3: Comparison of the local slopes of log(c) vs. log(t) between our model with 
p = 1/2 (*), and the anisotropic annihilation model, A + B — > (o). 
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Table 1 



p 


a 




7 


6 


V 





0.266(3) 


0.292(3) 


0.32(1) 


0.662(6) 


0.41(1) 


0.25 


0.289(3) 


0.311(3) 


0.34(1) 


0.660(7) 


0.39(1) 


0.5 


0.320(3) 


0.345(3) 


0.38(1) 


0.619(5) 


0.309(6) 


0.75 


0.380(2) 


0.400(4) 


0.44(1) 


0.585(5) 


0.215(7) 


1 


0.499(1) 


0.504(6) 


0.500(6) 


0.499(6) 


0.003(4) 



Table 1: Characteristic exponents for various values of the persistence parameter, p. The 
exponents are: c ~ t~ a , (Iaa) ~ (Iab) ~ t J , (L) ~ t s , (n) ~ t v . The error margins 
represent the spread of values obtained from extrapolations with various powers of 1/t. 
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